import matplotlib.pyplot as plt
import seaborn as sns
import pandas as pd
import sys
import numpy as np
import os
import loompy
import matplotlib as mpl
mpl.style.use('default')#classic
import scanpy as sc
from IPython import sys_info
print(sys_info())
seed = 123
np.random.seed(seed)
import scvelo as scv
np.random.seed(seed)
scv.logging.settings.logfile = "/myroot/03_merged/01_KNN/scv_log.txt"
scv.logging.print_version()
scv.settings.set_figure_params('scvelo')
scv.logging.settings.verbosity=2
exec(open('/myroot/03_merged/01_KNN/parameters.py').read())
print("n_neighbors_fate = "+str(n_neighbors_fate))
print("moments_neighbors = "+str(moments_neighbors))
print("moments_n_pcs = "+str(moments_n_pcs))
print("n_top_genes = "+str(n_top_genes))
loom files have been combined for all samples using:
loompy.combine(loom_files, output_file= output_file, key="Accession")
adata = scv.read(output_file)
Clusters information added from:
scv.pl.scatter(adata, color=['clusters'], alpha=.5,legend_loc='right margin',frameon=True)#
scv.pl.scatter(adata, alpha=.5,frameon=True,legend_loc='on data')
Basic preprocessing: filtering genes and normalizing data
scv.pp.filter_and_normalize(adata, min_counts=min_counts, min_counts_u=min_counts, n_top_genes=n_top_genes)
print("Parameters: "+str(min_counts)+str(min_counts)+str(n_top_genes))
Compute the first- and second-order moments (needed for velocity estimation)
scv.pp.moments(adata, n_pcs=moments_n_pcs, n_neighbors=moments_neighbors)
Estimation of velocities
scv.tl.velocity(adata)
Calculating probabilities of one cell transitioning into another cell, using cosine correlation and are stored in 'velocity graph':
scv.tl.velocity_graph(adata,random_neighbors_at_max=1000000)
scv.pl.velocity_embedding_grid(adata,basis='umap',scale=1,arrow_size=1,legend_fontsize=6,legend_loc='on data',color="clusters",arrow_length=1,density=1,frameon=True)
scv.pl.velocity_embedding_stream(adata, basis='umap',legend_loc='on data', show=True,alpha=.05,color='clusters')#
PAGA graph abstraction (Wolf et al 2019) provides a graph-like map of the data topology with weighted edges corresponding to the connectivity between two clusters. In scvelo, PAGA is extended by velocity-inferred directionality.
scv.tl.paga(adata, groups='clusters',use_rna_velocity=False)
scv.tl.paga(adata, groups='clusters',use_rna_velocity=True)
In new versions simply use the following to get and visualize the transition confidence information.
df = scv.get_df(adata, 'paga/transitions_confidence', precision=2).T
df.style.background_gradient(cmap='Blues').format('{:.2g}')
print(adata.uns["paga"]["transitions_confidence"][:2,:3])
Use diffrent transitions_confidence threshold for visualizations:
scv.pl.paga_compare(adata, basis='umap', threshold=0.3, arrowsize=30, edge_width_scale=.5,
transitions='transitions_confidence', frameon=True)
scv.pl.paga_compare(adata, basis='umap', threshold=0.2, arrowsize=30, edge_width_scale=.5,
transitions='transitions_confidence', frameon=True)
scv.tl.cell_fate(adata, n_neighbors=n_neighbors_fate,copy=False)
scv.tl.terminal_states(adata,copy=False,self_transitions=True)#,self_transitions=True)# default : self_transitions=False,
scv.pl.velocity_embedding_grid(adata,basis='umap',scale=1,arrow_size=1,legend_fontsize=6,legend_loc='on data',color="clusters",arrow_length=1,density=1,frameon=True)
scv.pl.velocity_embedding_grid(adata, color=[ 'root_cells', 'end_points'], legend_loc='on data',frameon=True)